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Abstract. Recently, secondary island formation due to the tearing instability of the 
Sweet-Parker current sheet was identified as a possible mechanism that can lead to fast 
reconnection (less sensitive dependence on Lundquist number S ) both in numerical 
simulations using Particle-in-Cell (PIC) method [Daughton et al. 2009], as well as us- 
ing resistive magnetohydrodynamics (MHD) [Lapenta 2008; Bhattacharjee et al. 2009]. 
This instability is thought to appear when S is greater than a certain threshold. These 
recent results prompt us to perform more resistive MHD simulations of a basic recon- 
nection configuration based on the island coalescence instability, using much higher 
resolutions and larger S . Our simulations are based on a fairly standard pseudo spectral 
code, which has been tested for accuracy, convergence, and compared well with codes 
using other methods [Ng et al. 2008]. In our simulations, formation of plasmoids were 
not found, except when insufficient resolution was used, or when a small amount of 
noise was added externally. The reconnection rate is found to follow the Sweet-Parker 
scaling when no noise is added, but increases to a level independent of S with noise, 
when plasmoids form. Latest results with S up to 2 x 10 5 will be presented. 



1. Introduction 

It is well known within the space physics community that the classic Sweet-Parker the- 
ory of magnetic reconnection predicts a reconnection rate that is too slow for most ener- 
getic phenomena, such as solar flares or magnetospheric storms. Strictly speaking, the 
Sweet-Parker theory is based on steady state two dimensional (2D) incompressible re- 
sistive magnetohydrodynamics (MHD), with a uniform classical resistivity tj, although 
it might also apply approximately in a more general system. It predicts a reconnection 
rate that is proportional to S~ 1 ^ 2 , where S = iiqVaLIt] is the Lundquist number, with Va 
being the Alfven speed based on the upstream magnetic field strength, and L being the 
typical length scale of system. Obviously such reconnection rate would be too small 
for large S , which is usually the case since L is very large for most systems in space or 
astrophysics, and n is usually small for high temperature plasmas. Such a small recon- 
nection rate is related to the formation of a long (to the order of the system size) but very 
thin (with thickness 8 also proportional to S _1 ' 2 ) reconnection current sheet, and so the 
inflow velocity is restricted due to the constrain of the continuity of mass. It is also well 
known that there have been many different theories or simulations over many years 
that can have a much higher reconnection rate, based on relaxing one or more restric- 
tions in the Sweet-Parker model, e.g., effects due to slow shock in compressible MHD, 
anomalous or non-uniform resistivity, Hall or two-fluid MHD, turbulence, and kinetic 
physics using direct particle-in-cell (PIC) simulations (too many to cite in this short pa- 
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per). In the last few years however, there has been a renewal of interest in reconnection 
in the Sweet-Parker setting, since a much higher reconnection r ate (apparently inde- 
pendent of 5) has been r e ported in some studies (lLapentall2008l: lLoureiro et al.1 2009 ; 
Bhattacharjee et al.|[2009l : ICassak et al.ll2009l : Huang & BhattacharjeeJ|2010j). A major 



reason for such an increase in reconnection rate is considered to be due to the tear- 
ing instability of the Sweet-Parker current sheet in the high-Lundquist number regime, 
and the formation and ejection of sec ondary islands (plasm oids), which has also been 
observed recently in PIC sim ulations (iDaughton et al.ll2009h . but actually has been ob- 



served for many years, e.g., (Lee & Ful ll986). Incidentally, one of the authors of this 
paper (Ng) has seen secondary island formation in reconnection simulations since 1996, 
but always regarded that to be due to insufficient resolution and thus such an observation 
was not reported. In this study, we plan to revisit this problem by performing careful 
simulations using as high resolution as we can (up to 8 192 2 ), for S up to 2 x 10 5 . The 
simulation geome try is as described in dNget al.ll2008l) and is essentially the same as 
what was used in (IBhattacharjee et al.ll2009h . What we found is that we have not been 
able to see the tearing instability of the Sweet-Parker current sheet, as well as the for- 
mation of plasmoids, desp ite the 5 value used is in the range that these are seen in 
(Bhattacharjee et al.ll2009l) . The reconnection rate found over a range of S follows the 
Sweet-Parker scaling of S~ l/1 . 

Moreover, in order to make sense of why we obtain different results, we performed 
another set of simulations with the same settings, except, random noise was now added 
in small scales. With that, plasmoids do appear and the reconnection rate found is 
roughly independent of S . This idea of reconnection with noise is similar to the concept 
of turbulent reconnection, which has been around for a while dMatthaeus & Lamkin 
1986; lLazarian & Vishniaclll999). There have a lso been a number of new studies on 
turbulent reconnection, e.g., ( Kowal et al.l I2009T). Note that some studies men tioned 
above also have noise added (lLoureiro et al.ll2009l: iHuang & Bhattacharieel l2010h . al- 
though there is still no universal understanding on the importance of the added noise. 

In Section |2j we describe briefly our simulation method and geometry. Our nu- 
merical results will be presented in Section |3] with a discussion in the final Section. 



2. Simulation Model 

Our results are based on numerical simulations of the 2D incompressible MHD equa- 
tions, which can be written in the following normalized form, 

d t Q. + [cf>, Q] = [A, J] + vV 2 ± Q, d,A + [0,A] = nV 2 ± A, (1) 

where A is the flux function so that the magnetic field is expressed as =z + V ± A x z, 
<p is the stream function so that the velocity field is u = V ± (p xz,Q. = -V 2 is the 
z-component of the vorticity, J = -V 2 ± A is the z-component of the current density, and 
[(p, A] = d y (pd x A - dyAdxCp. In this normalized form, the resistivity n is essentially the 
inverse of the Lundquist number. A more precise conversion of S = 2nA/n should be 
used when comparing with other studies, where A = 0.4 is the maximum magnitude of 
A used in the initial condition A(x,y, t = 0) = A sin(27r.x) sin(27ry). In all our runs, we 
have chosen viscosity v to be the same as 77. Periodic boundary conditions are used in 
both the x- and y- directions, with < x < 1 and < y < 1. A small flow is added 
initially to start the growth of the island coalescence instability. 
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The simulation code is based on the pseudo spectral method using fast Fourier 
transform (FFT) on a 2D bi-periodic domain. It is de-aliased by the standard 2/3 rule. 
The nonlinear term is calculated in the physical space on a uniform grid of colloca- 
tion points. A second order predictor-corrector method is used for time integration. 
The code is parallelized using a parallel version of the FFT. This code has be en tested 
for co nvergence, with accuracy of results compared with two other methods dNg et al.1 
2008) and found to be highly accurate in tests of conservation laws (with fractional 
errors usually three to four orders of magnitude smaller than those using a finite dif- 
ference method), as expected from a spectral code. However, the price to pay for a 
more accurate code is that it requires much longer CPU time to run, especially in high 
resolution cases due to the use of a uniform grid, as well as explicit time integration. 
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Figure 1 . (a) Left: linear resolution of adequately resolved runs for different re- 
sistivity j]. Dotted lines that scale with rj~ l l 2 and a dashed line that scales with 77 1 
are plotted for comparison, (b) Right: maximum reconnection rate over the whole 
duration of simulation for cases with different 77. Dotted lines scales with r] ' 2 . 



3. Simulation Results 

We have run many cases for different 77 and linear resolution (grid is at resolution 
square). Fig. 1(a) shows all cases (indicated with diamonds) that we have run and were 
adequately resolved. Note that for 77 at larger values, different resolutions have been 
used up to 2048 2 or 4096 2 so as to test for convergence. It is confirmed that as long as 
a run is adequately resolved, the numerical solution only change slightly as resolution 
increases. Moreover, we have also run cases with resolution lower by a factor of two of 
the lowest resolution cases shown in the figures. It usually was quite obvious that those 
runs were not well-resolved. Some of them indeed have plasmoids, which disappear as 
resolution is increased. Another remark on Fig. 1(a) is that the linear resolution seems 
to change from scaling with if ' 2 (as shown in dotted lines) in the large rj regime to 
scaling with if 1 (as shown in the dashed line) in the small 77 regime, with the tran- 
sition at around 77 ~ 10~ 4 . More careful examination of simulation outputs indicates 
that this is because at smaller 77, the length scale (proportional to 77) of other structures 
becomes smaller than the width of the Sweet-Parker current sheet, which scales with 
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T] l l 2 . Such small scale structures are found to be at the locations where the reconnection 
outflows collide with stationary structures. Moreover, these small scale structures have 
orientation perpendicular to the Sweet-Parker current sheet. Therefore, resolving the 
Sweet-Parker current sheet itself (say by using nonuniform grids with denser grid lines 
along the Sweet-Parker current sheet) is not enough for resolving the whole simulation 
domain. Obviously such a change of dependence of resolution on n means it is even 
more difficult to simulate high-Lundquist number cases, and one needs to pay attention 
to such a resolution demand when simulating these cases. 

In Fig. 1(b), we plot the maximum reconnection rate of each of these runs over the 
whole duration of simulation. Since we are interested in the high Lundquist number 
regime, we measure the reconnection rate simply by looking at the value of A at the 
central X-point, Aq, so that reconnection rate is given by dAo/dt. From this plot, it 
is clear that the data follow a Sweet-Parker scaling of n 1 ^ 2 , as indicated by the dotted 
lines. This means that we have not found increased reconnection rate even for the case 
with the smallest resistivity n = 1.25 x 10~ 5 , or S ~ 2 x 10 5 . More careful examination 
of the outputs, by plotting contours of A and J in movies, shows that there were no 
plasmoids forming in these runs. Fig. 2(a) shows Aq as functions of time t for cases 
with different S showing explicitly how Aq grows slower with t as S gets larger. 




Figure 2. (a) Left: Aq as functions of t for different S . (b) Right: Ay as functions 
of t for the S = 2.5 x 10 4 case using 1024 2 resolution with noise added (red) and 
without noise (black). 

Since no plasmoids were observed with the well-resolved runs, and that our ex- 
perience shows they can form in under-resolved cases, it is plausible that the inherent 
noise present in other numerical methods could have led to the formation of plasmoids 
in those simulations. To explore this possibility further, we have run our cases again 
with external random noise added throughout the whole simulations. The method o f 
adding noise is similar to the procedure outlined in ( Alveliuslll999l : lKowal et al] [2009). 



Noise is added to the simulation in spectral space with a characteristic wavenumber (de- 
pending on the system parameters) and an amplitude that ensures that the noise added 
is always of a small fraction of the kinetic and magnetic energies at any time. 

Fig. 2(b) shows an example of this study with Aq vs t plotted for the S = 2.5 x 10 4 
case using 1024 2 resolution with noise added (red) and without noise (black). Obvi- 
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ously one can see that the reconnection rate is enhanced in the case with noise added. 
An examination of outputs shows plasmoid formation and ejection, e.g. as shown in 
Fig. [3j so that instantaneous reconnection rate fluctuates a lot at the same time, but 
averaged out to be a higher rate. 
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Figure 3. (a) Top: Contour lines of A for the case S = 1.25 x 10 4 with no noise. 
Color filled contours are for J to show the location of the current sheet (dark blue 
region), (b) Bottom: plot for the same case with noise added when a plasmoid is 
formed. The two plots are at roughly the same time. 




Figure 4. (a) Left: Aq as functions of t for different S with noise added, (b) Right: 
average reconnection rate dAo/dt for different S with noise added (red diamonds) 
and without noise (black squares). 

In Fig. 4(a), Aq as functions of t for more cases with different S are plotted on 
the same graph, showing that the reconnection rate does not change much even as S 
changes a lot. These cases are simulated using different resolutions, and thus the scale 
at which noise is added also changes. However, to have a more consistent comparison, 
such noise is kept to have energy roughly of the same ratio to the total energy (~ 10~ 6 ), 
and that the correlation time of the random noise is kept at about 0.25 for all cases. To 
characterize reconnection rates for all these cases, the averaged values of dAo jdt over 
the connecting phase are plotted on Fig. 4(b) in red diamonds, together with reconnec- 
tion rate for the original cases without added noise (black squares). This shows clearly 
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that the cases with noise have a reconnection rate roughly independent of S , and are at 
a level much larger than the reconnection rate without noise, especially in the large S 
limit, since the reconnection rate is consistent with a Sweet-Parker scaling of S~ l/2 . 

4. Discussion and Conclusion 

We have presented our simulation results showing that reconnection rate in resistive 
MHD do follow the Sweet-Parker scaling over a significant range of Lundquist number 
up to S ~ 2x 10 5 , if runs are well-resolved, and that no external noise is added. In these 
runs, apparently no tearing instability occurs at the Sweet-Parker current sheet, and thus 
no plasmoid formation is observed. In trying to understand why our results differ from 
results from other studies, we repeated our runs with random noise added externally. 
We have shown that plasmoids do form in these cases, resulting in a reconnection rate 
roughly independent of S, and thus much higher than the case with no noise in the 
large S limit. While these results do show that plasmoid formation in reconnection of 
this setting occurs only when external noise is added, we do not claim that this also 
explains what were observed in simulations of other studies. This result does highlight 
the importance of making sure that a simulation is run with sufficient resolution, e.g., 
by checking convergence with higher resolutions. 

We remark that we have not done a full study of how reconnection rate depends 
on properties of the added noise, which is in fact not an easy task. For example, one 
surely can expect that the effect of noise would be negligible if its total energy tends to 
zero. We have confirmed this by running cases at different levels of noise energy, but 
we have not run enough cases so far to show a dependence of reconnection rate on the 
noise energy. More importantly, the next step we need to investigate is what physical 
effects are responsible for generating such noise that can enhance reconnection rate. 
Can we stay within MHD, e.g., using turbulence, or do we need some external physics? 
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